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We construct a family of many-body wave functions to study the many-body localization phase 
transition. The wave functions have a Rokhsar-Kivelson form, in which the weight for the config¬ 
urations are chosen from the Gibbs weights of a classical spin glass model, known as the Random 
Energy Model, multiplied by a random sign structure to represent a highly excited state. These 
wave functions show a phase transition into an MBL phase. In addition, we see three regimes of 
entanglement scaling with subsystem size: scaling with entanglement corresponding to an infinite 
temperature thermal phase, constant scaling, and a sub-extensive scaling between these limits. Near 
the phase transition point, the fluctuations of the Renyi entropies are non-Gaussian. We find that 
Renyi entropies with different Renyi index transition into the MBL phase at different points and 
have different scaling behavior, suggesting a multifractal behavior. 


I. INTRODUCTION 

For an isolated quantum many-body system totally de¬ 
coupled from the environment, the system can act as its 
own heat bath. One way to distinguish systems which 
reach the thermalized state from those which don’t is by 
studying the entanglement entropy of a subsystem. If a 
subsystem is in thermal equilibrium, the entanglement 
and thermal entropies must be the same and thus must 
satisfy a volume law, i.e. , the entropy should scale like 
Ni in d dimensions. This is different from the typical 
scaling behavior of the entanglement entropy of ground 
states, which under some constraints, satisfy an area law 
scaling.liH^ 

Evidence has recently accumulated for a certain class 
of interacting systems with quenched disorder which fail 
to thermalize. These systems go under the name many- 
body localizatiorPH^ (MBL) and are the interacting ana¬ 
logue of Anderson insulators.!^ This failure to thermal¬ 
ize has been attributed to an extensive number of lo¬ 
cally conserved charges.lSHIllFor a review of MBL phases, 
see the recent review of Ref. M In an Anderson in¬ 
sulator, all single-particle eigenstates are exponentially 
localized in real space, quantum diffusion is impossible 
at zero temperature and the system is an insulator on 
macroscopic scales.!^ In such a state thermalization is not 
possible (without an external heat bath). It has been sug¬ 
gested that many-body localization re sults from a similar 
localization of states in Hilbert space.^^E^ 

The phenomenon of MBL is principally (and theoret¬ 
ically) observed in states deep in the excited state spec¬ 
trum of a macroscopic system, and hence have an ex¬ 
tensive excitation energy which, following standard (but 
somewhat inexact) terminology, we will call ‘finite energy 
density’ states. As a function of some tuning parameter 
(typically disorder), there can be a phase transition from 
an ergodic (thermalized) phase to an MBL phase. In 
contrast to its non-interacting counterpart (the Ander¬ 
son insulator), an MBL phase transition can occur at 
finite temperature.!^ We should note that phonons may 


interfere with the observation of an MBL phase in solids, 
but MBL may be physically realized in optical lattice 
systems, see e.g. Ref. [161 

In this paper we consider the problem of the MBL 
phase transition by constructing an ensemble of simple 
‘model’ many-body wave functions with a simple struc¬ 
ture parameterized by a ‘disorder’ strength, and study 
the phase transition as a function of this parameter. We 
will show that, in spite of their simple structure, these 
model wave function can represent both thermal states 
and MBL states. The wave functions that we consider 
have a structure similar to the Rokhs ar-Kivelson (RK) 
stateJi^ and their generalizationsMore specifically, 
we consider states that are linear superpositions of quan¬ 
tum states labeled by the classical configurations of a 
system of N Ising spins, and have the form 

\^rem) = Y.W[C]\C), (1.1) 

{C} 


Here the quantum mechanical amplitude W[C] for a con¬ 
figuration C of the Ising spins is given by the Gibbs weight 
of a classical spin glass model known as the Random En¬ 
ergy Model, i.e. 

W[C] cx (1.2) 


where the ‘energy’ E[C\ assigned to the configuration C is 
taken to be a random number drawn from a Gaussian dis¬ 
tribution. By construction, the amplitudes of these states 
W[C] are positive real numbers. The associated classical 
spin glass model in infinite space dimension (since each 
spin is coupled to all the other A^ — 1 spins) is known 
to have a classical thermodynamic phase transition to a 
spin glass state.^^The parameter /3, which in the classical 
spin glass model is the inverse temperature, but in this 
work will be used as a parameter of the wave function. 
Notice that we have not defined a quantum Hamiltonian 
for which the wave function of Eq.(l.l| is an eigenstate 
and, hence, we have not actually defined an energy for 
the quantum system. Thus, the ‘energy’ of the Random 
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Energy Model should not be confused with the energy of 
the quantum state. 

The quantum state I'i'R EM) has the manifestly posi¬ 
tive weights shown in Eq.(1.2). Such a state can be a 
natural candidate for a ground state of a Hamiltonian 
but not for a typical excited state whose amplitudes are 
generally non-positive. To mimic a ‘typical’ state deep 
in the spectrum of a quantum system, we generalize this 
construction so that, for a given configuration C, the am¬ 
plitudes for these new states are just a random sign mul¬ 
tiplied by the amplitude W{C) discussed above (a similar 
approach has been used in Ref HH and Refill We will 
denote the new wave functions by |'I'/j£;M-i-sign[/3])- Here 
we will also consider the wave function without random 
signs denoted by and compare the physical 

properties of both types of wave functions. 

An advantageous aspect of our approach is that we 
have more analytical control over this system then is 
typical in interacting disorder systems. In addition we 
are able to perform numerical calculations with a system 
size (> 30) which can only be achieved in other num erical 
MBL studies with the use of matrix-product state^^HSl!, 
These wave functions are also conceptually simple, mak¬ 
ing them an ideal setting to further our understanding of 
MBL. 


Laumann, Pal and Scardicchio studied numerically the 
MBL state in the quantum Random Energy Model and 
found that the MBL quantum phase transition is dis¬ 
tinct from the quantum phase transition to the spin-glass 
phase.^ Here we will find that in the RK wave functions 
l^flBM-i-sign), although they are not actual eigenstates of 
the quantum REM model, the MBL and spin-glass tran¬ 
sitions are also separate. 

In this work we will be focused on three particular as¬ 
pects of the problem- the transition to the MBL phase, 
the scaling of the entanglement entropy with subsystem 
size and the transition from being geometrically delo¬ 
calized to localized in the Hilbert space. To identify the 
ergodic and MBL phases, we use the Renyi entanglement 
entropies Sn (where n is the Renyi index ) 


Sn = -^\ogTvp\ (1.3) 

1 — n 


of a subsystem A whose size Na < N/2 is smaller than 
half of the entire system. In the limit n —>■ 1, Sn con¬ 
verges to the von Neumann entropy. In an ergodic sys¬ 
tem, there is a regime where the Renyi entropy obeys a 
volume law which changes linearly as a function of Na- 
Bauer and Nayak have argued that, for most states in an 
MBL phase, their entanglement entro]w scales at most 
as an area law of the subsystem size.l^ We will show 
below that, in our MBL phase, the Renyi entropy is sub- 
extensive as a function of subsystem size and is bounded 
by a finite constant deep inside the MBL phase (Fig. 
[^. We take particular note of volume laws at an energy 
density that corresponds to infinite temperature (ITV) 
which scales as Na log 2. The scaling behavior of the en¬ 
tanglement entropy may depend on the subsystem size 


often showing a crossover from ITV to sub-extensive as 
the subsystem gets larger. An important subtlety of our 
model comes from the lack of correlation length in the 
REM (inherent in a system at infinite dimension). This 
makes the relevant parameter to consider in looking for 
such a crossover not the absolute size Na but the ratio 
t = Na/N and our results will be quoted as a function 
of this parameter. 

While entanglement entropy can be used to distin¬ 
guishing MBL from ergodic phases, geometric localiza¬ 
tion is a measure of compactness of the wave function in 
Hilbert space. In an Anderson insulator, the localization 
is of the single-particle wave function in real space and 
can be characterized by the inverse participation ratio 
defined as 


Yn= I d‘'x|V'(x)|2" 


(1.4) 


where \'ip{x)\'^ is the probability distribution of single¬ 
particle state in real space.^^ Generically, Yn takes the 
scaling form 


Yn - (1.5) 

where the exponent T(n) = Dn{n — 1). For the exten¬ 
sive (delocalized) state, Dn = d, while for the localized 
state Dn = 0. For the critical single-particle wave func¬ 
tion at the mobility edge, Dn has a non-trivial depen¬ 
dence on n, and it indicates that the critical wave func¬ 
tion has a multifractal nature.l^S^^ ^ For non-interacting 
systems, this multifractal behavioilS^I^ is also manifest 
in the Renyi entropy for the single-particle critical wave 
function.lsami 

Here we will present evidence that multifractal be¬ 
havior also appears in the entanglement near the phase 
transition into the MBL phase by looking at the scal¬ 
ing behavior of the Renyi entropies. Loosely speaking 
the many-body generalization of this multifractality mea¬ 
sures the degree of localization of states in the multi¬ 
dimensional Hilbert space (in a real-space basis) and not 
just those of a single particle orbital. Multifractal be¬ 
havior of weight of a state in a Hil bert space has been 
studied recently by several authors these stud¬ 
ies multifractality is used to characterize the geometry of 
a state in Hilbert space, i.e. its degree of localization in 
the Hilbert space. In those works, the Shannon-Renyi en¬ 
tropies used to quantify the multifractal behavior of the 
many-body states is a measure of the statistical proper¬ 
ties of the wave functions as probability distributions and 
are unrelated to the concept of quantum entanglement. 
Multifractality has also been discussed in connection with 
the fidelity of the ground state wave functions in systems 
at the infinite-disorder fixed point.l^ So far as we know, 
multifractal behavior of quantum entanglement in wave 
functions close to the MBL transition has not been dis¬ 
cussed previously in the literature. This is one of the 
main questions that we address in our work. 
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A. Summary of the main results 


(a) REM model 

Non-Gaussian fluctuation 


In thi s paper, w e consider the MBL transition, the 
value of S{Na,N) = {S{Na, N)) and the geometric lo¬ 
calization of the wave function |'I'i?BM-i-sign[/3])- We show 
that the wave-function geometrically localizes in Hilbert 
space at (3 = 1.18. W e also find three regimes of en¬ 
tanglement scaling for S2{Na',N = oo ): a thermalized 
regime (at infinite temperature) where the entanglement 
entropies for a partition show volume law scaling as a 
function of the subsystem size Na] a regime bounded 
by constant entanglement entropy; and a regime which 
is sub-extensive but not constant. The transition from 
sub-extensive to constant happens at, or before, the ge¬ 
ometric localization transition. Notice that the presence 
of extensive scaling for any Na > 0 implies that the finite 
size scaling of S 2 {N;Na) as a function of TV at a hxed 
Na is also extensive. Because S 2 {pa) is hard to compute 
exactly, we analytically calculate lower and upper bounds 
for it. In addition, we identify the MBL transition (with 
respect to the second Renyi entr opy) v ia a numeric al scal - 
ing collapse with (S' 2 ) and ((5S'2]P^^ (see Section IVA). 
We find that the the MBL transition is different from the 
geometric localization transition and therefore the wave- 
function is still geomet rically delocalized in Hilbert space 
at the MBL transitiorl^^ l^ . While the entanglement 
entropy S 2 {Na,N = 00 ) scales sub-extensively with Na 
at the MBL transition, the transition to sub-extensive 
scaling doesn’t correspond to the MBL transition. In 
addition, the numerical evidence suggests that the MBL 
transition (for S' 2 , at Na/N = 1/2) happens at a /3 where 
S 2 {N-, Na/N = 1/2) as a function of N still scales exten¬ 
sively. While these statements hold for S 2 , we also con¬ 
sider Sn, finding that the MBL transition as well as the 
analytic bounds for sub-ITV scaling happen for different 
n at different p. The former of these appears to scale 
linearly in (n — l)/n. In addition, the scaling exponents 
identifled from scaling collapse are different for different 
n. This is an indication of multifractal behavior. We fur¬ 
ther study the entanglement spectrum in this regime and 
observe an entanglement gap between the lower continu¬ 
ous band and the other higher eigenvalues and show that 
this feature explains transitions which differ for different 
n. In Figjl] a and b, we give a succinct and broad sum¬ 
mary of the phase diagram (including the phase diagram 
of the classical REM for comparison). 

We also explain the importance of the random sign 
structure. This effect we compare the behavior of the 
second Renyi entropies {S 2 {pa)) for |'Lij_EM-i-sign) with 
I'^rem) without random sign. The difference, AS' 2 , is 
found to decrease monotonically as a function of disor¬ 
der strength. In the low disorder regimes, (S' 2 (pa)) for 
I'^rem) is a constant and thus their difference AS 2 takes 
the maximal value. At the values of t < 1/3, the differ¬ 
ence disappears before entering into the MBL phase. 

The rest of the paper is organized as follows. In Section 
|TTj we explain the construction of a Rokhsar-Kivelson- 
type wave function which assigns to the amplitude of a 
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FIG. 1. (Color online), (a) Phase diagram for the classical 
REM model as a function of j3. The spin glass phase tr ansition 
occurs at /3 = \/2 log 2. The shaded region between \/log 2/2 
and yj2 log 2 has non-Gaussian fluctuation, (b) Schematic 
phase diagram for the REM wave function in terms of the 
entanglement entropy showing different scaling behaviors in 
four different regimes as a function of /3. In the regime (i) 
{Sn) is equal to Na log 2 for all Na- In the intermediate 
regime (ii), {S 2 } is sub-extensive but does not saturate to a 
constant. Regime (iii) is the MBL phase. In this regime, the 
wavefunction is not localized in the Hilbert space, (iv) is also 
the MBL phase with the wavefunction localized in Hilbert 
space. 


quantum many-body state the Gibbs weight of a classi¬ 
cal spin glass model (the Random Energy Model) and 
introduce a random sign structure into it to access repre¬ 
sent typical excited states with a finite excitation energy 
density. In Section [III A| we briefly review the classical 
random energy model (REM) and interpret the classi¬ 


cal spin glass phase transition in it. In the Section HIB 


we construct a REM wave function with a random sign 
structure to mimic a highly excited state for a Hamil¬ 
tonian with disorder. We will analytically compute the 
Renyi entropy to show that there is a thermalized regime 
and a many-body localized phase. In Section|IV]we calcu¬ 
late the Renyi entropy numerically and find the location 
of MBL phase transition by finite size scaling. In the 
Section |IVD| we study the Renyi entropy for the REM 
wave function without random sign and demonstrate the 
importance of the sign structure. In the Section |V] we 
summarize our results and conclude that there is a MBL 
phase transition in the REM wave function with random 
sign structure. 


II. ROKHSAR-KIVELSON MODEL WAVE 
FUNCTIONS 


In the physics of strongly correlated systems there are 
many examples in which the properties of a new state of 
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matter can be represented by simple ‘model’ wave func¬ 
tions. The best known examples of such model wave 
functions include the BCS wave function for the ground 
state of a superconductor and the Laughlin wave function 
for the fractional quantum Hall fluid. 


model, and the quantum critical point in the wave func¬ 
tion corresponds to the classical phase transition at tem¬ 
perature l//3c- 


A. RK-type wave functions 

In this paper we will give a description of the MBL 
states using a particularly simple class of model wave 
functions with a structure similar to the one proposed 
by Rokhsar and Kivelson (RK) to capture the physics 
of the ground states of strongly frustrated quantum 
antiferromagnets.^In the RK construction the quantum 
mechanical amplitude of a many-body state is given by a 
local function of the degrees of freedom, as expressed in 
the orthonormal basis set {|C)} where C is the ‘configura¬ 
tion space’. In the RK problem the configuration space 
is the set of dimer coverings of a 2D lattice. While in 
the RK case the dimers are a qualitative representation 
of spin singlets on each bond of the lattice, a picture of 
this type has been generalized to many other systems, in¬ 
cluding Kitaev’s Toric Code state,!^ which represents the 
topological (or deconfined) phase of a Z 2 gauge theory. 

Since the weights of the RK wave functions are lo¬ 
cal and positive, they can also be regarded as the Gibbs 
weights of a related problem in classical statistical me¬ 
chanics with the same degrees of freedom on the same 
lattice. Thus if the basis of orthonormal states is the 
set {|C)}, i.e. such that {C\C') = Sc.c, the generalized 
normalized RK states are 

= ( 2 . 1 ) 

IC} 


where 


Z = (2.2) 

{C} 

Here E[C] and Z are, respectively, the energy for the 
classical configuration C and the partition function for 
the associated classical problem, and /3 plays the role of 
the inverse temperature. The original RK wave function 
is an equal-amplitude superposition of all dimer cover¬ 
ings in 2D. This wave function can be associated with 
the partition function for the classical dimer model and 
is also the ground state of the quantum dimer model at 
the critical poiniP^tHl (for 2D bipartite lattices) and a 
Z 2 topological stat^^ (for non-bipartite lattices). States 
of these type are exact ground states of a special type of 
quantum Hamiltonians that are the sum of projectio n op- 
erators and are closely related to classical dynamics. 

The generalized RK wave function \'^rem) inherits 
many properties from the classical model. For instance, 
the equal time correlation function of the wave function 
is the same as the correlation function of the classical 


B. Random sign wave function 


We want to consider our family of wave functions as 
representing finite energy density states with the possi¬ 
bility of supporting both ergodic and MBL phases. As 
we noted above, there exists Hamiltonians constructed 
by projection operators for which the RK wave functions 
of the form of Eq.( 2.1) are the exact ground states. How¬ 
ever, we should note that its amplitudes are all strictly 
positive. If the classical Hamiltonian used to generate 
the RK state only has local interactions, we can show this 
state must have an area law. This follows from the fact 
that one can write a Schmidt decomposition of the state 
where the number of terms is bounded by the number of 
classical interactions which are broken (See Appendix A 
for detail). The entanglement properties of these type of 
states have been discussed in weat detail in the case of 
the quantum dimer model,^^®^ of the associated quan- 
tum Lifshitz model,and of the Toric Code statel^^*^ 


To have the potential for ergodic states, then, we must 
either use a non-local classical Hamiltonian or introduce 
directly a more rich sign structure into the wave function 
as done, e.g., in Refs. m and While the classical 
Hamiltonian we are using for the REM model is non¬ 
local, it nonetheless supports an area law at low disorder; 
in fact, at /3 = 0, the entanglement entropy is zero over 
any cut. Therefore, we introduce a random sign structure 
in the wave function giving 


\'i>REM+sign) =sc e (2.3) 

{C} 


where sc is a quenched random sign associated to each 
configuration C. The average ratio between the number 
of positive and negative signs is one. A similar construc¬ 
tion was discussed in Ref. [551 where they showed that the 
random sign structure can lead to a thermalized phase. 
In this paper, we will study the MBL phase transition in 
a RK-wave function with random sign. 


III. MANY BODY LOCALIZATION PHASE 
TRANSITION 


In this section we will show that the RK wave function 
with random signs has an ergodic (thermalized) regime 
and an MBL phase. We begin with a summary of the 
properties of the classical Random Energy Model whose 
Boltzmann weights will enter into the structure our wave 
function. 
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A. The Random Energy Model 

The random energy model (REM) is a simple classi¬ 
cal model which has a phase transition to a spin glass 
phasePSl This effectively infinite dimensional model has 
2^ configurations and is the infinite range coupling limit 
of the Sherrington-Kirkpatrick spin glass modell^In the 
REM model, the energy for each configuration is no 
longer given by the complicated spin glass Hamiltonian, 
but rather simply an independent random variable. This 
random variable has the Gaussian distribution 

P{E) = {2N-K)-^/'^e-^. (3.1) 

The number of configurations in the energy interval 
[fVe, 7V(e -I- ^)], in expectation, is 

(AA(e, (e + ,5))) = (3.2) 

where 

E 

s{x)=\og2-—, witha;=— (3.3) 

In the thermodynamic limit, iV —^ oo, this expectation 
value takes the asymptotic form 

lim (A/'(e, (e -k 5))) = exp{fV max^,^;^ e+ 5 ]s(a;)} (3.4) 

N^oo 

Notice that s(x) > 0 only in the interval x G [—e*,e*], 
where e* = \/2 log 2. This means that for e S [—e*,e*], 
(A/’(e)) is exponentially large, and the fluctuations are 
very small. For e outside the interval [—e,,e*], A/’(e) is 
exponentially small. 

The partition function for this model is simply given 

by 


2 

Z = = f dEM{E) = f dx 

(3.5) 


where ^(x) equals to 


(/)(x) = log 2 - y - /3x, (3.6) 

Similar to the calculation for (A/’(e)), we can also use 
the saddle point approximation to calculate the partition 
function to obtain 


Z = exp{A max[(()(x)]} (3.7) 


which is the exact result in the thermodynamic limit, 
iV —>■ oo. By computing ^i)„iax = max[(/)(x)], it is easy to 
show that the free energy density equals to 


/(/?) 


log Z _ (/fniax 

PN ~ p 

_P _ ]og2 R R 

2 P F ^ Psg 

-(2 log 2 ) 1/2 


(3.8) 


At Psg = •\/21og2, there is a discontinuity in the second 
derivative of the free energy density, showing that there is 
a phase transition at this point, the spin glass transition. 

This phase transition can be further stu died b y com¬ 
puting the inverse participation ratio (IPR) ,1^2111] 


2 y2 nPE. 

y^{P)^Y.p" = 




(3.9) 


This quantity measures how many configurations effec¬ 
tively contribute to the partition function and measur¬ 
able quantities. The expectation value for Y 2 equals tcP 


{y2m 


0 P<Psg 

1 ~ ^ y ^ Psg 


(3.10) 


At low temperatures, P > Psg, the participation ra¬ 
tio takes a finite value. This means that the system 
is completely frozen to 0 ( 1 ) number of configurations 
and is in a non-ergodic phase. At high temperature, 
P < Psg, all configurations will contribute to the ther¬ 
modynamic properties of the system. For instance, in 
the high temperature limit /3 —> 0, the Boltzmann’s mea¬ 
sure becomes uniform and the second participation ratio 
becomes Y 2 {P —>■ 0) = 2“^. When iV —^ 00 , F 2 —>■ 0. In 
general, when P < psg, Yn scales with the system size 


~ V = 2^ (3.11) 

From the results on logZ(/3), we can compute the expo¬ 
nent T(n) (for n > 1) 

({n-l){l-jn), 0 < 7 <;^ 

r(n) = < n(l - ^)2, ^<7<1 (3-12) 

[ 0 , 7>1 


a2 

where 7 = ^ upper bound 

for the Renyi entropy of the wave function |ff'/j£;M+sign) 
later. 

Although the REM model is a simple toy model, it 
has a spin glass phase transition, i.e., it undergoes a lo¬ 
calization transition. It also shows a rich structure in the 
fluctuation of the free energy. According to the results of 
Ref. [SH for the ergodic phase, which occurs for P < Psg, 
there are two regimes with different fluctuation behavior 
of the free energy. When P < y^log2/2, the fluctuations 
of the free energy are Gaussian, and satisfy the central 
limit theorem. When Psg > p > ^J\og2/2, there are 
non-Gaussian fluctuations of the free energy driven by 
the Poisson process of the extreme values of the random 
energies. This regime does not satisfy the central limit 
theorem. The resulting phase diagram of classical REM 
model is shown in Fig{^ (a). 
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B. Random sign REM wave function 


2. Analytic Bounds 


We can now construct a quantum state following the 
procedure of Eq.(|2.3|), 


\^rem) = (3.13) 

fc} 


We now analytically compute bounds to the Renyi en¬ 
tropy. Notice that for a disordered system, we need to 
take a quenched ensemble average of the Renyi entropy, 

{Sn{pA)) = p^j,) (3.15) 

1 — n 


whose amplitudes are the Boltzmann weights of the clas¬ 
sical REM. By further introducing the random sign struc¬ 
ture, the wave function takes the form 

I'i’REM+sign) = (3.14) 

fC} 

Here both sq and E[C\ are random variables and are in¬ 
dependent of each other. The random sign sc takes the 
values ±1 with equal probability over the entire Hilbert 
space of 2^ spin configurations. 

Before we do any calculations, we can first estimate the 
scaling behavior of the Renyi entropy for this quantum 
state in the extreme limits. In this wave function, /3 is 
a tuning parameter and has the physical meaning of the 
disorder strength. The random sign structure is used to 
represent a highly excited quantum state. When there is 
no disorder, i.e., /3 = 0, the amplitude for every configu¬ 
ration is the same. At this point, the Renyi entropy for 
the subsystem A with different Renyi index is equal to 
the thermal entropy at infinite temperature. Thus this 
wave function is thermalized, and the entanglement en¬ 
tropies obey a volume law.^^ As /3 increases, the disorder 
becomes stronger and the entanglement entropy becomes 
smaller. The wave function is eventually localized to a 
small fraction of the configurations as /3 > Psg- As we 
already showed in the previous section, since the number 
of these configurations is only 0(1), in this regime the 
Renyi entropy is bounded by a finite constant. 


1. Nomenclature and thermodynamic limit for scaling 

The entanglement entropy (Sn) at fixed (3 depends on 
both the subsystem size Na as well as the total system 
size N and there can be separate functional forms for the 
scaling of {Sn) with respect to either of these parameters. 
While in many physical systems, these scalings coincide, 
this is not the case in our model and so it is important to 
be clear about the distinction. Unless otherwise specified, 
our discussion will always focus on scaling with Na at 
fixed system size N. 

For the thermodynamic limit considered in this work, 
we let both Na and N go to infinity but their ratio t = 
Na/N to be a finite value. The changes in entanglement 
scaling we identify then happen at particular values of t 
for 0 < t < 1/2 for different values of /?. The MBL phase 
in our terminology is identified with scaling collapse. 


In this quenched average the disorder is frozen and does 
not evolve with time. We mainly focus here on the second 
Renyi entropy, but the results are easily extended to the 
other Renyi entropies. 


We begin by noting that the quenched average in 
Eq.(3.15) is redundant if, in the thermodynamic limit, 
the system is self-averaging.l^ However, this will not al¬ 
ways be the case. For brevity, we sometimes will denote 
(^„(p^)) = (Sn). 


Let us define a bipartition our system of N spins into 
two subsets (or regions), A and B. The reduced density 
matrix for region A is 


P 


A 

a,a' 



\Ca) ® {Ca'\ 

(3.16) 


where p"^ is the unnormalized reduced density matrix. 

For the above reduced density matrix, {Sn{pA)) equals 
to 


(SniPA)) = 


1 — n 


((logTrpl)-n(logZ(/3))) (3.17) 


The second term (log Z(/?)) can be calculated by the sad¬ 
dle point approximation and the result is already shown 
in Eq.(3.8). However, the first term (logTrp)!) is hard to 
obtain analytically. Instead of calculating it directly, we 
compute a lower bound and upper bound for it. 

Lower Bound : Consider the annealed average 


SuHpa)) = (log(Tr/^) - n(logZ(/3))) (3.18) 

which is much easier to compute. By using Jensen’s 
inequality,^ when n > 1, it is straightforward to see that 
the annealed average of the Renyi entropy Sn{{pA)) pro¬ 
vides a lower bound for the quenched average {Snip a)), 

i.e. 


{SniPA)) > SniipA)) (3.19) 

To obtain the annealed average S 2 i{pA)), we need to 
calculate (Trp^). By means of simple manipulations we 
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find 


(Tr^A)=E((?t')^) 

a,a' 

= E((^a,a)')+E((^a.aE) 

a a^a' 

a b 

+ E^(E®“’'>®“' 

a^a' b 

+ (2^B _ 1)2^6"^^' + (2"^^ - 1)2^6"^^' 

(3.20) 

where the last step is derived by using that for a Gaussian 
distribution (e““^) = e“ 

In the thermodynamic limit, when Na < N, the an¬ 
nealed average of the second Renyi entropy becomes 


-§(Ea.,b+E^ 


S2i{pA)) = 


fiVAlog2, /3 </3i 

\7V(log2-/32), I3>p, 


(3.21) 


where /3i = y^(l — t) log 2 and t is the ratio Na/N. The 
result is plotted in Figj^ (a), where the black dashed 
curve is the lower bound for {S 2 )/St at t = 1/3. St = 
Na log 2 is the thermal entropy for subsystem Na at in¬ 
finite temperature. Notice that when /? < /3i, since the 
lower bound S' 2 ((aa)) = St, {S 2 [pa)) must be equal to 
fV^log2. When /? > v^log2, S 2 {{pa)) < 0, this lower 
bound is replaced by zero and is not useful anymore. 

Upper Bound : It is easy to see that Trp^ can be 
bounded from below by the second participation ratio of 
the classical REM. Indeed, 


TrpA > 




> EM =F2(/3)(3.22) 


Z (^)2 Z (/ 3)2 

Thus, an upper bound for the quenched average {S 2 {pa)) 
is given by 



FIG. 2. (Color online) (a) The lower and upper bound for 
(S 2 ). The black dashed curve is for the lower bound and 
the red solid curve is for the upper bound. The system size 
N = 300 and the ratio t = 1/3. (b) The upper bound for 
(S„). The setup is the same as (a). 


This upper bound indicates that when 


/3 >/3ub{n) = {1 - \ -— -t)^/2Ao^ (3.27) 

V n 


{S 2 {pa)) <T2iP)Nlog2 


(3.23) 


where T2(/3) is defined in Eq.( 3.12l . This upper bound 
puts a constraint on {S 2 {pa)), showing that when 

/? > /3„6(n = 2) = (1 - y|)v/21og2 (3.24) 


then 


(^2(P.4)) <^.4log2 (3.25) 


which is the thermal entropy at infinite temperature. 
The red solid curve in Eigj^ (a) is the upper bound for 
{S 2 )/St at t = 1/3. 

Using a similar approach, we can prove that when n 
is an even number, Trp)^ > Yn{f5), and the quenched 
averaged nth Renyi entropy satisfies 


{Sn{pA)) < 


r„(/3)Nlog2 
n — 1 


(3.26) 


we obtain the bound 

{Sn{pA)) < NA\og2 (3.28) 

which, again, is the thermal entropy at the infinite tem¬ 
perature. The result for the upper bound is shown in 
Eigj^ (b). Notice that when /3 < /3„b(n = 00 ), the differ¬ 
ent curves in Figm (b) are overlapping with each other. 
When /3uf,(n = ooj^< /3 < figg, the upper bounds for (S'„) 
show different scaling behavior. Combining the upper 
and lower bounds for different n leads to a regime where 
the n = 2 Renyi entropy satisfies a volume law and be¬ 
haves like a thermal entropy at infinite temperature, but 
for n > 2 is strictly below this bound. 

a. Implications from the bounds : From the above 
results, we find an upper and lower bound for the Renyi 
entropy. Here we consider (S' 2 ) in detail, which thus has 


































Intermediate regime (S2) 


same for both states, and it is given by the probability 
of the configuration C in the classical REM, 



P 


FIG. 3. (Color online) A summary of our knowledge of the 
phase diagram of REM+aign, based on (Sn). The dotted 
black line indicates a bound to the left of which {S2} analyt¬ 
ically follows the T = 00 volume law (ITV). The dotted red 
line indicates a bound to the right of which {S2) analytically 
is strictly below the T = 00 volume law. j3sg indicates an¬ 
alytically the transition to a localized Hilbert space for any 
Na/N and a guaranteed constant for all {Sn}- The solid pur¬ 
ple dots indicate the numerically computed transition points 
for {S2) (the purple line is a fit for the eye). In comparison, 
the blue dashed line indicates a bound to the right of which 
(Soo) analytically is strictly less then ITV, and suggests mul- 
tifractality. 


the following behaviors 

r {S2) = S 2 {{pa)) = Na log 2 , /3 < /3i 

I 7V(log 2 - /32) < {S 2 ) < Na log 2, P 2 < ^ < Psg 
[(^ 2 ) <-log(I-/3«g//3), /3>/3«g 

(3.29) 

where 132 = I3ub{n = 2) = (1 - log 2. 

For (5'2), there is a regime where (5'2) satisfies the 
volume law and equals the thermal entropy at T = 00 . 
When /3 > /3sg, (S' 2 ) is bounded by a finite constant. A 
phase transition into the MBL phase is expected to be 
between /?i < /3 < f3sg ■ To find the location of the MBL 
phase transition, we will use scaling collapse in Section 

CYl 


P[C] = \{C\'i>REM)r = -e 


2 = ±p-3ElC] 
Z 


(3.30) 


Given this probability distribution, we can compute its 
Shannon-Renyi entropies, which one can immediately 
see, c.f. Eq.(3.9), to be the same as the IPR of the classi¬ 


cal REM. Hence, the localization of wave function in the 
Hilbert space can be characterized by the IPR defined in 
the configuration space. For REM wave function. 


Yn = Y.mc)f^ = 

{c} (E 


— n^Ei 


^2^ 


-0Ei 


(3.31) 


This is the same as the IPR for the classical REM model 


defined in Eq.(3.9). Since T(n), defined in Eq.Q.ll I, is 


the multifractal spectrum, given explicitly in Eq.(3.12|, 


in the regime 0 < < Psg, the wave function itself has 

multifractal behavior. The multifractality of the wave 
function indicates the pre-freezing behavior before enter¬ 
ing into the MBL phase. Similar phenomenon was also 
observed in Refsl31|and|351 where they found that at the 
MBL phase transition point, the whole wave function is 
still delocalized in the configuration space. 

We can now use the results of the inverse participa¬ 
tion ratios of the classical REM summarized in section 


HI A to draw conclusions on the degree of localization 
in the 2^-dimensional Hilbert space of the I'^^rem) and 
EaBM-i-sign) wave functions. From the results of Section 


HI A we find that for /3 > /3sg = -^2 log 2 all the inverse 
participation ratios are finite as N ^ 00 and, hence, 
that the Shannon-Renyi entropies for the wave functions 
are finite (and are not extensive). Thus, in this regime 
these wave functions are (exponentially) localized. On 
the other hand, for f3 < Psg, the IPRs of the REM vanish 
exponentially fast as —>■ 00 , and so do the Shannon- 
Renyi entropies of the wave functions. In this regime the 
wave functions are not localized. In this regime the mul¬ 
tifractal nature of these wave functions is manifest in the 
size dependence of their Shannon-Renyi entropies. We 
should emphasize that, beyond setting the bounds that 
we have discussed earlier in this section, the knowledge 
the behavior of the Shannon-Renyi entropies alone yields 
no information on the scaling of quantum entanglement, 
which is our main interest. 


C. Localization properties of the wave function 

We end this Section with a discussion of the statistical 
properties of the states and |4'j?£;M+sign)- Given 

these states, we can define the amplitudes, i.e. their over¬ 
lap with an eigenstate of the spins \C). For each state, the 
square of the amplitude defines a probability distribution 
for the configuration C to occur in the state (and hence 
the random sign does not affect the probability distribu¬ 
tion). The resulting probability distribution is thus the 


IV. NUMERICAL RESULTS ON MBL PHASE 
TRANSITION 

In the previous section, we were able to analytically es¬ 
tablish the existence of a regime in which the quenched 
averaged second Renyi entropy {S 2 ) is strictly less than 
ITV scaling, and happens strictly before the localization 
transition. Nothing prevents this regime from being one 
in which {S 2 {Na)) still scales linearly but at finite T; 
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in addition, it doesn’t separate the MBL transition from 
the localization transition. To establish this, we numer¬ 
ically identify the transition point via scaling collapse. 
We explicitly construct different disorder samples from 
the random sign REM wave function and calculate the 
quenched average {Sn) with at least 1000 disorder realiza¬ 
tions. First, we focus on the n = 2 case, and then discuss 
the other values of n. As a sanity check, we verify that 
our numerical results are consistent with the analytical 
results above at large and small p. 

The following is a summary of our numerical results. 

At small P : Analytically we anticipate ITV for 
n = 2. In Fig. we concretely consider p — 0.3 at 
N = 30 with 2000 disorder realizations. We see that it 
corresponds to the expected (S' 2 ) = iVAlog 2 . 

As Na log 2 is the maximal entanglement entropy for 
any given realization, for a T = 00 volume law to hold, 
essentially all but a measure zero fraction of configura¬ 
tions must have this entropy. From the inset of Fig|^, 
one can see that the standard deviation of {S 2 ) at /3 = 0.3 
is zero showing this is indeed the case. Similar behavior 
can be observed for (Sn) with other Renyi indices, where 
{Sn) = Na log 2 and (SSn) is close to zero as /? < /?*. 

At large P > P^g '■ We find analytically that (5'„) 
with n > 2 is bounded by a finite constant and, hence, 
in this range the state is in the MBL phase. Note that a 
constant for small n gives a bound for the entanglement 
entropy for all larger m as Sn > Sm if to > n. In Fig|^ 
the numerical results for /? = 1.5 are presented, and show 
that when Na increases, (Sn) saturates to some constant 
value. Notice also from the inset of Fig[^, that as we 
move deeper into the localized phase, the standard devi¬ 
ation is monotonically decreasing. 

The change from ITV to constant in the entanglement 
entropy can be seen in Figj^ which shows d{S 2 )/dNA 
as a function of p. As system sizes increases, the slope 
quickly approaches log 2 for P < 0.5, and approaches 
0 for /3 > 1.2. The slopes for all 5 values of N start 
to dro p ar ound P = 0 . 6 , which is less than the P 2 = 
(1 — a/1/6)-v/ 2 log 2 « 0.6967 in Eq.(3.29) for subsystem 
ratio t = 1/3 and is consistent with the analytical result. 


A. Finite-size scaling 

To locate the MBL phase tra nsitio n, we use both {S 2 ) 
and its standard deviation, ^ 5*2 ) 42 | 43 | instance, in the 
bottom inset of Fig. it is shown that when P > 0.5, 
5 S 2 /ST will increase rapidly and reach the maximum 
value at some p. Fig. (bottom) shows that the stan¬ 
dard deviation at intermediate values of P is actually a 
non-trivial fraction of the maximum thermodynamic en¬ 
tropy. In fact, because of the breakdown in the REM of 
the Central Limit Theorem, the fluctuation of the Renyi 
entropy is non-Gaussian. This is seen, for example, in 
the distribution of S 2 /ST at P = 0.72 shown in Figj^ 
which is peaked around 1 and has a long tail. This long 
tail has power law scaling behavior shown in the inset 
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FIG. 4. (Color online) The numerical results for (Sn) with n 
from 1 to 00 at different values of /3. The total system size is 
N = 30. Each point in the plot is averaged over 2000 disorder 
realizations. At /3 = 0.3, S„ shows ITV behavior for all n. At 
P = 0.66,1.1, (Sn) deviates from the ITV with increasing /3, 
but the speed at which it moves away depends on its Renyi 
index. At P = 1.5, (Sn) is a constant for all n. 

of Fig0 with a power-law exponent between -3 and - 
2 , which implies a well-defined average but an infinite 
variance in the thermodynamic limit (which is consistent 
with the peak scaling as St)- This may be related to the 
quantum Griffiths phase found in Ref. saiMi and 1661 



P 


FIG. 5. (Color online) vs P graph at t = 1/3, where the 
slopes are obtained from finite differences. As P increases, 
the slopes drop from log 2 towards 0 for all 5 different system 
sizes. 

While we can use the peak of 5 S 2 and the transition of 
{S 2 ) from maximal to zero to locate the phase transition, 
these quantities scale with system size. We therefore use 
finite-size scaling to find the MBL phase transition. In 
the regime of interest, we perform a scaling collapse of the 
data for the two ratios {S 2 )/St and 6 S 2 /ST separately. 
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0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 1 

82/5^ 


(a) 



Xn=(P-Pc)N“ 


(b) 

FIG. 6. (Color online), a) Scaling collapse of {S2) / St at 
t — 1/3, where St = A^Alog2. It can be noticed that b 
is very close to 0. b) Scaling collapse of 5S2/St- The left 
insets show the original curves. Two scaling collapses give 
very close /3c values. Error analysis is only performed on the 
scaling collapse of {S 2 )/St- 


using scaling functions of the form 

7V'’$((^ -/3c)iV“) (4.1) 

and determine the form of the scaling functions <I>(ai) nu¬ 
merically. The exponent b is expected to be very close to 
0 for {S 2 )/St- When doing the scaling analysis, we find 
that the quality of the {S 2 ) /St collapse is better than 
that of the 582 /St- Error analysis is only performed on 
the scaling parameters obtained from {S 2 )/St, although 
generally the 882 /St collapse yields similar values for /3c. 

Notice that this scaling form assumes the existence of 
only one transition in spite of the fact that we have ana¬ 
lytical bounds that show the presence of ITV, constant, 
and sub-ITV scaling; numerically the latter appears to 


FIG. 7. (Golor online) The distribution of S2/St at /3 = 0.72. 
The total number of samples is 2000. The number of bins used 
to make the histogram is 100. The inset shows a log-log plot 
with a linear fit, which indicates a power law distribution. 
The number of bins in the log-log scale is 20 for S2 / St from 
0 to 1. However, as many of the bins have 0 sample, some 
points are not included the log-log plot. 


be neither constant nor linear with subsystem size. 

Fig. ID shows the scaling collapsed for (S' 2 ) at t = 1/3. 
In previous work on MBL phase transition, u = a~^ 
in Eq.( |4.1[ ), is the critical exponent for the localiza¬ 
tion length and is expected to satisfy Harris inequality 
1 / > 2/d, where d is the spatial dimension. However 
REM model is a highly non-local model and d here is 
equal to infinity, which suggests that there is no bound 
for ly. To make sure the scaling parameters obtained are 
truly associated with this transition, one can use them 
to scale the (S'„)/S't data following Ref. 011 For each 
/3, we plot (Sn) vs N as shown in the inset of Fig. , 
from which one can see that the scaled curves bifurcate 
smoothly into the two branches depending on /3. Fur¬ 
thermore, Fig. ID indicate that the system experiences a 
phase transition rather than a crossover, because of the 
clear separation of the two branches at large N. 

We have also done scaling collapse for (S' 2 ) with other 
subsystem ratios. At t = 1/2, we get /3c = 0.693; and 
at t = 1/5, we get /3c = 0.965. These results are shown 
as purple dots on Fig. [D All our collapsing results are 
summarized in Table IH 


Based on the above numerical results, we conclude that 
the MBL transition (both of (S 2 ) and all the other (S„) 
(see Sec |IV C and Appendix C) happen separately from 
the localization transition at /3sg- The MBL phase tran¬ 
sition is at Pc and is smaller than Psg. Between Pc and 
Psg, although it is in the MBL phase, the whole wave 
function is still delocalized in the Hilbert space. Only 
when P > Pgg, the wave function becomes localized in 
the Hilbert space and is already deep in the MBL phase. 
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FIG. 8. (Color online) This graph uses the scaling parameters 
obtained from Fig. Imto scale the {S 2 )/St vs N curves for /3 
values ranging from 0 to 2. The scaled graph clearly shows two 
branches - (1) curves with large j3 flow to low entanglement 
entropy; (2) curves with small /3 flow to T = oo entanglement 
entropy, i.e. NA^og2. The right inset shows the original 
curves, where each curve corresponds to a different f3. 



t 

n 

a 

b 

P. 

{Sr.) 

S'jp 

1/3 

1 

0.811 

0.007 

0.970 

1/3 

1.5 

0.897 

0.0093 

0.8680 

1/3 

2 

0.906 

0.0006 

0.8072 

1/3 

3 

0.902 

0.0095 

0.755 

1/3 

20 

0.8389 

0.025 

0.6657 

1/3 

CK) 

0.817 

0.039 

0.651 

1/2 

2 

0.90 

0.070 

0.693 

1/5 

2 

0.766 

0.010 

0.965 

5Sn 

St 

1/3 

1 

0.893 

0.321 

0.983 

1/3 

1.5 

0.902 

0.396 

0.861 

1/3 

2 

0.923 

0.400 

0.798 

1/3 

3 

0.908 

0.399 

0.753 

1/3 

20 

0.835 

0.383 

0.655 

1/3 

CK) 

0.824 

0.397 

0.643 

1/2 

2 

0.909 

0.090 

0.678 

1/5 

2 

0.856 

0.147 

0.963 


TABLE I. Collected scaling collapse data 


B. Intermediate regime 


According to the scaling collapse in Fig. the MBL 
phase transitio n for {S 2 ) happens at around 0.8 at t = 
1/3. From Eq.(3.24), we also know that when /3 > 0.697, 
{S 2 ) < Na log 2. This implies that there is an intermedi¬ 
ate regime between ITV and MBL phase (regime (ii) in 
Figj^(b)). In this regime, (S' 2 ) is sub-extensive and the 



xn=(P-P,)n‘' 


(a) 



XN=(P-Pe)N“ 


(b) 

FIG. 9. (Color online) (a) Scaling collapse of {Scx,)/St at 
t = 1/3, where St = 77a log 2. It can be noticed that b is 
very close to 0. (b) Scaling collapse of SSoa/Sr- The left 

insets show the original curves. Two scaling collapses give 
very close Pc values. Error analysis is only performed on the 
scaling collapse of {Soo)/St- 


slope 


0 < s{Na) = 


1 d{S2) 

log 2 PNa 


< 1 . 


(4.2) 


This regime is not a cross-over but sharply defined with 
non-analyticities in the curve s{t) signaling its beginning. 
While a combination of our analytical bounds and numer¬ 
ical results can bound the location of this transition, we 
are not able to numerically pinpoint its location. It is 
interesting to note, though, that all the curves in Fig. 
seem to cross at a single point; for t = 1/3, this point 
is approximately /3 = 0.72 which is surprisingly close to 
the P 2 bound. While we can’t say anything dehnitive 
about the value of the s{Na) in the intermediate regime, 




















































C. Multifractality of the Renyi entanglement 
entropies 
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FIG. 10. (Color online) This graph uses the scaling parame¬ 
ters obtained from Fig. |^to scale the (Soo) / St vs N curves 
for P values ranging from 0 to 2. The scaled graph clearly 
shows two branches - (1) curves with large /? flow to low en¬ 
tanglement entropy; (2) curves with small /3 flow to T = oo 
entanglement entropy, i.e. Na log 2. The right inset shows the 
original curves, where each curve corresponds to a different 
/?• 


t=l /3 



While we have identified transitions in (S' 2 ), we can 
also consider (Sn) for n ^ 2. The analytical bounds 
(Figj^ show that there is a regime where S 2 is still ITV 
where S'oo scales at a rate less then ITV. We can also use 
scaling collapse to identify the MBL phase transition in 
(Sn)', for example, see Fig. and 10 More scaling col¬ 
lapse graphs with different n and subsystem ratio t can 
be found in Appendix C. We summarize the N^/N = 1/3 
results for various n in Fig |ll[ from which one can clearly 
see the transition depends on n. We attribute this as a 
sign of multifractal behavior similar to the multifractal 
behavior found in the critical wave function of the Ander¬ 
son localization phase transition point .1^ In the Anderson 
localization problem it is known that multifractality is a 
feature of the wave function for the mobility edge.^It is 
not known if this is also the case in MBL or if there is a 
multifractal phase. From our data we cannot at present 
make a definitive determination. It is interesting to note 
that j3c{n) is linearly proportional to (n—l)/n. The phase 
transition point /3c for von Neumann entropy is 0.97 and 
we identify this as the true MBL transition. 


We can understand the different scaling behavior of 
{Sn) by considering the entanglement spectrum. Fig(^ is 
the distribution of eigenvalue A of pA for a randomly cho¬ 
sen disorder configuration at different (3. When /3 = 0.3, 
A forms a continuous band around 0.001 which is ap¬ 
proximately equal to 1/2^-^. While for /3 = 0.66,1.1 and 
1.5, there is an obvious gap between the lower continuous 
band and the other higher eigenvalues. The entanglement 
gap increases as /3 increases. 


This inspires us to write down a simplified two-level 
model for entanglement spectrum which includes a flat 
band and a single Xmax- If we assume that (Sn) is ap¬ 
proximately equal to its upper bound in Eq.( |3.26 1 , by 
using Sn=oo = - ^Og Xmax, we have 


Xmax — ‘2 


-( 1 - 


✓2 log 2 


pN 


(4.3) 


FIG. 11. (Golor online) /3c vs curve at t = 1/3. The 
linear fit line’s slope is —0.3213 ± 0.0034, and its y-intercept 
is 0.9709 ±0.0023. 

Fig.j^shows only two plateaus suggesting that the s{Na) 
in this intermediate regime is not constant. Instead we 
conjecture that the slope changes continuously as a func¬ 
tion of Na in this regime in a nonlinear way. This means 
that this regime is non-thermal and doesn’t correspond 
to a thermal density matrix at any temperature. It is 
not clear whether to call this regime a separate phase, 
particularly as the data collapse on {S 2 )/St and 6 S 2 /ST 
only identify a single transition. 


All the other Ai « 1/2^'* and there is a gap between Ai 
and Xmax- This toy model exhibits the multifractal be¬ 
havior described above with a transition to constant en¬ 
tanglement slope that scales with n > 2. The transition 
in this simplified model is direct from ITV to constant; 
to capture our intermediate regime, the single Xmax can 
be replaced by a finite number Xi with each of them are 
separated by a finite gap. Also the lower flat band can be 
replaced a continuous band with more complicated band 
structure. These additional ingredients are required to 
have an intermediate regime with sub ITV as well as 
accurately finding the constant n = 1 von Neumann en¬ 
tropy. 
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N=30, t=l/3, 1000 samples 



FIG. 12. (Color online) Entanglement spectrum for one ran¬ 
domly chosen disorder configuration at /3 = 0.3, 0.66,1.1 and 
1.5, with N = 30, Na = 10. The gaps in the entanglement 
spectrum become evident when /3 gets larger. 


D. The random sign structure in the wave function 


We have introduced a random sign structure to convert 
our ground state wave function into one at finite energy 
density. While for any strictly positive wave function, 
the introduction of random signs can only increase the 
entanglement entropy, it is interesting to ask what effect, 
if any, the random sign has here. This requires com¬ 
puting the Renyi entropy {Snip'^)) for the REM wave 
function without the random sign |'I'_r_em)- Different 
from I'l'flBM-i-sign), the Renyi entropy for \^rem) is not 
a monotonic function with j3 and has more complicated 
scaling behavior. In fact rem ) has zero entanglement 
entropy at both /3 = 0 and /3 —>■ oo. It is interesting 
to note that while {S 2 {pa)) and {S 2 {p'a)) have long tails 
and infinite variance in the thermodynamic limit at in¬ 
termediate (3, their difference AS 2 = {S 2 {pa)) — {S2{p'a)) 
appears to have finite variance (as shown in Fig. |13[ ). 

Since the quenched average ( 5*2 (^( 4 )) is hard to access 
analytically, we can compute a lower bound via Eq.(3.18) 
using Jensen’s inequality^ (see Appendix B for details). 
We find that if t < 1/3, there is a region log 2 < (3 < 
— log 2 where \'^rem) has ITV entanglement. As 
this is the maximal allowed value, it then directly follows 
that there is no difference between {S 2 ) and (S' 2 ) due to 
the introduction of signs. 

Moreover, we conjecture that AS '2 decreases monoton- 
ically as a function of [3. This is consistent with the nu¬ 
merical results shown in Fig |14| Following from this con¬ 
jecture, we would have that AS 2 = 0 for all (3 > y/2t log 2 
and t < 1/3. This is because, for all t < 1/3 both models 


show ITV between f3 = y'2t log 2 and a/(1 — t) log 2 and 
hence AS 2 = 0. The regions B, C, D in Fig. [T5| denote 
where \'^ rem ) and |'I'_R_EM-i-sign) have the same (S' 2 ). 


FIG. 13. (Golor online) {S 2 ) of REM-fisign and REM, and 
their difference. The 1000 disorder configurations are the 
same for REM-(-sign and REM. Error bars here represent 
standard deviations. The difference of the (S 2 ) has well 
bounded variance for all p. 



P 


FIG. 14. (Color online) Difference of {S„) between REM-|-sign 
and REM for n = 0.1,...,0.9 in 0.1 steps (triangles), and for 
n = 1.5, n = 2,..., 20, and 00 (full circles), at A = 30 and 
Na = 10. The number of disorder configurations is 2000. The 
random configurations are not the same for REM-|-sign and 
REM. 


Having identified regimes where the introduction of 
random signs doesn’t affect the entanglement entropy, 
we also identify regimes where the entanglement entropy 
can be shown to be different. When (3 = 0, I'^'rem) is a 
constant (actually a product state) whereas I'k/ESM+sign) 
is a volume law. We can argue that this extends to larger 
p. Defining X = Tr/i^i and Y = Tr(p( 4 )^, where pA is 
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S 2 (REM+sign) lower bound ^ 



No random signs, N=30, (1=0.4 



FIG. 15. (Color online) Different behaviors of AS 2 = 
{S 2 {REM + sign)) — {S 2 {REM)) on the t vs j3 graph. Region 
A represents the area where {S 2 {REM + sign)) follows ITV 
and {S 2 {REM)) is a constant, so A5'2 obeys ITV. Region B 
represents the area where {S 2 {REM + sign)) and {S 2 (REM)) 
are both ITV, and AS 2 is zero. As a result, A5'2 is also zero 
in Region C based on the monotonicity argument. In Region 
D, both {S 2 {REM + sign)) and {S 2 {REM)) are localized in 
Hilbert space and AS '2 is zero. We do not have enough infor¬ 
mation to determine the behavior AS '2 in Region E. 


the unnormalized reduced density matrix, we have 


(log^)<log(^)<log^ (4.4) 


where the second inequality, while not true in general, 
appears to be numerical validated in our case. In the 
thermodynamic limit, log can be directly computed 
(see Eq.(3.20) and Eq.(Bl)). We find that when t < 1/3 
and /3 < ^/2t log 2, log =0. In this region, (S 2 ) for 
l^flBM-i-sign) continues growing as volume while rem) 
stays constant, indicating that the random sign structure 
can thermalize the wave function and is responsible for 
the volume law scaling behavior. Fig. 16 is the numerical 
result for (S/) at /3 = 0.4. We can see that when n > 2, 
they all saturate to a constant. This region is highlighted 
in blue on Fig. |15| and marks a region where the models 
differ maximally. 


Finally, we note that ASn can be numerically com¬ 
puted. Fig{^ is ASn for iV = 30 and Na = 10. We 
find that for all n, ASn = 0 when (3 > 0.8. This re¬ 
sult, while only for iV = 30 and so not absent finite-size 
effects happens to be at the location of the {S 2 ) MBL 
phase transition. 


FIG. 16. (Golor online) The numerical results for {S'„) with¬ 
out random sign with n from 1.5 to 00 at /3 = 0.4. The total 
system size is N = 30. Each point in the plot is the average 
of 2000 disorder realizations. 


V. CONCLUSIONS 


We have studied the many-body localization phase 
transition in a class of many-body wave functions. We fo¬ 
cused our analysis in a class of wave functions, \^rem), 
whose amplitudes are the Boltzmann weights of a clas¬ 
sical spin glass model in infinite dimension, the Random 
Energy Model. In order to mimic the structure of wave 
functions of highly excited states with a finite excita¬ 
tion energy density we considered another class of states, 
l^'flBM-i-signs), whose amplitudes are obtained by multi¬ 
plying the amplitudes of \^rem) by a random sign for 
each configuration. We studied the MBL problem in the 
|^'fi£;M-i-signs) wave function, for different regime of the 
parameter /3, by using both analytical and numerical ap¬ 
proaches. 

We showed that there is a direct phase transition into 
the MBL phase. Here we assume that the MBL phase 
is characterized that the entanglement entropies, as a 
function of the size of the observed region Na, scale to 
a constant value, a feature that we observed explicitly 
for large enough values of j5. The location of the phase 
transition point is identified by scaling collapse of Renyi 
entropy and its standard deviation. In the thermalized 
regime, there is a regime where the Renyi entropies with 
different Renyi index all equal to the thermal entropy 
at T = 00 . When (3 > Pc, the system enters into 
the MBL phase where the entanglement entropy is sub- 
extensive. The MBL phase transition point Pc is smaller 
than Psg and this suggests that the MBL phase transi¬ 
tion and the classical spin glass phase transition are dif¬ 
ferent. Upon entering the MBL phase, the random sign 
structure is not important any more and {ASn) between 
l^'fiBM-i-sign) and I'i’REM) is zero. For /? > Psg, the wave 
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function is deep inside the MBL phase and only 0(1) 
number of the configurations in the REM wave function 
contributes significantly to the statistical average. The 
Renyi entropy in this regime will reduce to a finite con¬ 
stant. We find that close to the phase transition point 
/3c, the fluctuation of Renyi entropy is strong in the fi¬ 
nite size system. In this regime, Renyi entropies with 
different Renyi index show different scaling behavior and 
are similar to the multifractal behavior observed at the 
Anderson localization phase transition point. Moreover, 
(Sn) has a phase transition at different /3c. Finally we 
note that we have refrained ourselves from performing 
the same extensive studies for the wave function without 
random signs, While we have evidence that this 

wave function too has a thermalized regime, since it has 
only strictly positive amplitudes we do not expect it to 
provide an useful description of the MBL problem. Nev¬ 
ertheless it may be useful to investigate its properties in 
a separate publication. 

Note added in the proof: After this work was accepted 
for publication we became aware of the recent indepen¬ 
dent work of Yang, Chamon, Hamma and Mucciolc!^ 
who have also discussed the spectrum of the reduced den¬ 
sity matrix in several models, including the wave func¬ 
tions discussed in this work. 
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Appendix A: Entanglement entropy for the 
RK-wave function with classical local Hamiltonian 


For the RK-wave function defined in Eq.(2.1), if the 
related classical model has a local Hamiltonian, the en¬ 
tanglement entropy satisfies the area law. This scaling 
behavior only relies on the property of the RK state and 
does not depend on whether the quantum model is crit¬ 
ical or not. Here, we will briefly review the calculation 
following RefUBJ 


For the RK state with a classical local Hamiltonian, 
after partitioning the system into two parts A and B, it 
can be approximately written in this way, 

|4/o)=^Ar|4'^(r))|4/B(r)) (Al) 

r 


where 41^ (P) and ^'^(r) are the RK wave functions de¬ 
fined in region A and region B with the same boundary 
configuration P 


l^'A(r)) = E 

.4 

B 


VzW) 

^/ZbW) 


\ca(T)} 

\cBiT)) 


(A2) 


When r ^ P', they satisfy (4>A(r)|'kA(r')) = 0 and 
('I'B(r)l'I'B(r')) = 0. The summation in Eq.(Al) is the 
sum over all possible boundary configurations along the 
cut and Ar = ^/Zb^/VZ. 

Thus Eq. ( |A1[ ) is the Schmidt decomposition of the 
wave function, the reduced density matrix in regime A is 


PA = Y.xl\'SA{r)){^A{r)\ (A3) 

r 


Since the dimension of pA only depends on the dimension 
of the Hilbert space along the boundary, the entangle¬ 
ment entropy should satisfy the area law. 


Appendix B: Lower bound for {S2{p'a)) 


For the REM wave function without random sign, the 
quenched average {S 2 {p'a)} is hard to access analytically, 
instead we calcul ate the annealed average S 2 {{p'a}) de¬ 
fined in Eq.( 3.18| , which gives the lower bound for (S”^). 
To obtain S^IXpa))^ need to know {Trfp'^)'^) hrst, 
where p'^ is the unnormalized reduced density matrix. 


(Tr(p'A)2) = 2 ^ 6 "^^^^ + (2^® - 1)2^6^^^" 

-k(2^^ - I)2^e^^' -k (2^-^ - I)(2^® - l)2^e^ 

(Bl) 

This gives 5*2 ((p^)) in the thermodynamic limit. When 
0 < t < 1/3, there are three regimes, 


/3<y2nop 

'S' 2 ((pa)) = N A^Alog2, V2ilog2 < <)log2 

[N{\og2-P^), /3 > V(I - <)log2 

(B2) 

When 1/3 < t < 1/2, there are two regimes 


S2i{p'A)) 


N{\og2-/3^), 


/3< V21og2/3 
/3> V21og2/3 


(B3) 



















16 


^ L. Bombelli, R. K. Koul, J. Lee, and R. D. Sorkin, Phys. 

Rev. D 34, 373 (1986). 

^ M. Srednicki, Phys. Rev. E 50, 888 (1994). 

® J. Eisert, M. Cramer, and M. B. Plenio, Rev. Mod. Phys. 
82, 277 (2010). 

^ D. M. Basko, I. L. Aleiner, and B. L. Altshuler, Ann. Phys. 
321, 1126 (2006). 

® V. Oganesyan and D. A. Huse, Phys. Rev. B 75, 155111 
(2007). 

® A. Pal and D. A. Huse, Phys.Rev.B 82, 174411 (2010). 

^ P. W. Anderson, Phys. Rev. 109, 1492 (1958). 

® M. Serbyn, Z. Papic, and D. A. Abanin, Phys. Rev. Lett. 
110, 260601 (2013). 

® D. A. Huse and V. Oganesyan, Phys. Rev. B 90, 174202 
(2014). 

J. Z. Imbrie, “On Many-Body Localization for Quantum 
Spin Chains,” (2014), arXiv:1403.7837 

V. Ros, M. Muller, and A. Scardicchio, Nuclear Physics 
B 891, 420 (2015) 

A. Chandran, 1. H. Kim, G. Vidal, and D. A. Abanin, 
Phys. Rev. B 91, 085425 (2015) 

R. Nandkishore and D. A. Huse, Annual Review of Con¬ 
densed Matter Physics 6, 15 (2015) 

C. Monthus and T. Garel, Phys. Rev. B 81, 134202 (2010). 
E. Canovi, D. Rossini, R. Fazio, G. E. Santoro, and 

A. Silva, Phys. Rev. B 83, 094431 (2011). 

W. Morong and B. DeMarco, Phys. Rev. A 92, 023625 
(2015) 

D. S. Rokhsar and S. A. Kivelson, Phys. Rev. Lett. 61, 
2376 (1988). 

E. Ardonne, P. Fendley, and E. Fradkin, Annals Phys. 
310, 493 (2004). 

E. Fradkin, Field Theories of Condensed Matter Systems 
(Cambridge University Press, Cambridge, UK, 2013). 

B. Derrida, Phys. Rev. Lett. 45, 79 (1980). 

T. Grover and M. P. A. Fisher, J. Stat. Mech. , PlOOlO 
(2014). 

T. Grover and M. P. A. Fisher, “Entanglement and 
the Sign Structure of Quantum States,” (2014), 
arXiv: 1412.3534, 

V. Khemani, F. Pollmann, and S. Sondhi, “Obtaining 
highly-excited eigenstates of many-body localized Hamil¬ 
tonians by the density matrix renormalization group,” 
(2015), arXiv: 1509.00483 

X. Yu, D. Pekker, and B. K. Clark, “Finding ma¬ 
trix product state representations of highly-excited eigen¬ 
states of many-body localized Hamiltonians,” (2015), 
arXiv: 1509.01244, 

C. R. Laumann, A. Pal, and A. Scardicchio, Phys. Rev. 
Lett. 113, 200405 (2014) 

B. Bauer and C. Nayak, Journal of Statistical Mechanics: 
Theory and Experiment 9, 5 (2013) 

D. J. Thouless, Physics Reports 13, 93 (1974). 

F. Evers and A. Mirlin, Rev. Mod. Phys. 80, 1355 (2008). 
®® F. Wegner, Zeitschrift fiir Physik B Condensed Matter 36, 

209 (1980). 

C. de C. Chamon, Mudry, and X. G. Wen, Phys. Rev. 
Lett 77, 4194 (1996). 

H. E. Castillo, C. C. Chamon, E. Fradkin, P. M. Goldbart, 
and C. Mudry, Phys. Rev. B 56, 10668 (1997). 

I. A. Gruzberg, A. W. W. Ludwig, A. D. Mirlin, and M. R. 


Zirnbauer, Phys. Rev. Lett 107, 086403 (2011). 

V. E. Kravtsov, 1. M. Khaymovich, E. Guevas, and 
M. Amini, “A random matrix model with two localization 
transitions,” (2015), arXiv:1508.01714 
T. C. Halsey, T. Jensen, L. Kadanoff, 1. Procaccia, and 

B. Shraiman, Phys. Rev. A 33, 114 (1986). 

M. H. Jensen, L. Kadanoff, and 1. Procaccia, Phys. Rev. 
A 36, 1409 (1987). 

X. Jia, A. R. Subramaniam, 1. A. Gruzberg, and 
S. Ghakravarty, Phys. Rev. B 77, 014208 (2008). 

S. Ghakravarty, Int. J. Mod. Phys. B. 24, 1823 (2010). 

X. Ghen, B. Hsu, T. L. Hughes, and E. Fradkin, Phys. 
Rev. B 86, 134201 (2012). 

D. J. Luitz, F. Alet, and N. Laflorencie, Phys. Rev. Lett. 
112, 057203 (2014) 

E. J. Torres-Herrera and L. F. Santos, Phys. Rev. B 92, 
014208 (2015) 

R. Vasseur and J. E. Moore, Phys. Rev. B 92, 054203 
(2015) 

J. A. Kjall, J. H. Bardarson, and F. Pollmann, Phys. Rev. 
Lett. 113, 107204 (2014). 

R. Vosk, D. A. Huse, and E. Altman, Phys. Rev. X 5, 
031032 (2014). 

A. D. Luca and A. Scardicchio, Europhys. Lett. 101, 37003 
(2013). 

D. J. Luitz, N. Laflorencie, and F. Alet, Phys. Rev. B 91, 
081103 (2015). 

A. Y. Kltaev, Annals of Physics 303, 2 (2003). 

R. Moessner and S. L. Sondhi, Phys. Rev. Lett. 86, 1881 

( 2001 ). 

C. Henley, Journal of Physics: Condensed Matter 16, S891 
(2004). 

C. Castelnovo, C. Chamon, C. Mudry, and P. Pujol, An¬ 
nals of Physics 318, 316 (2004). 

S. Furukawa and G. Misguich, Phys. Rev. B 75, 214407 
(2007). 

C. Castelnovo and C. Chamon, Phys. Rev. B 76, 184442 
(2007) 

J. M. Stephan, S. Furukawa, G. Misguich, and 
V. Pasquier, Phys. Rev. B 80, 184421 (2009). 

E. Fradkin and J. Moore, Phys. Rev. Lett. 97, 050404 
(2006). 

B. Hsu, M. Mulligan, E. Fradkin, and E.-A. Kim, Phys. 
Rev. B 79, 115421 (2009). 

B. Hsu and E. Fradkin, J. Stat. Mech. P09004 (2010). 

A. Hamma, R. lonicioiu, and P. Zanardi, Physics Letters 
A 337, 22 (2005). 

M. Levin and X.-G. Wen, Phys. Rev. Lett. 96, 110405 
(2006). 

V. Khemani, A. Chandran, H. Kim, and S. L. Sondhi, 
Phys. Rev. E 90, 052133 (2014). 

D. Sherrington and S. Kirkpatrick, Phys. Rev. Lett. 35, 
1792 (1975). 

B. Derrida, Phys. Rev. B 24, 2613 (1981). 

M. Mezard and A. Montanari, Information, Physics, and 
Computation (Oxford University press, Oxford, UK, 2009). 
A. Bovier, 1. Kurkova, and M. Lowe, Ann. Prob. 30, 605 
( 2002 ). 

E. Buffet, J. Phys. A 26, 1823 (1993). 

M. Reed and B. Simon, I: Functional Analysis, revised 
ed.. Methods of Modern Mathematical Physics (Academic 


33 

34 

35 

36 

37 

38 

39 

•41 

42 

43 

44 

45 

46 

47 

48 

49 

50 

51 

52 

53 

54 

55 

56 

57 

58 

59 

60 

61 

62 

63 

64 





17 


Press, London, UK, 1981). 

K. Agarwal, S. Gopalakrishnan, M. Knap, M. Mueller, and 
E. Demler, Phys. Rev. Lett. 114, 160401 (2015). 

A. C. Potter, R. Vasseur, and S. A. Parameswaran, Phys. 
Rev. X 5, 031033 (2015) 

Z.-C. Yang, C. Chamon, A. Hamma, and E. R. Mucciolo, 


“Two-component Structure in the Entanglement Spectrum 
of Highly Excited States,” (2015), arXiv: 1506.01714 
J. M. Stephan, G. Misguich, and V. Pasquier, JSTAT-J. 
Statist. Mech.: Theor. Exp. 2012, P02003 (2012). 




